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ABSTRACT 

Context. Understanding the growth of the cores of giant planets is a difficult problem. Recently, Lambrechts & Johansen 
(2012; LJ12) proposed a new model in which the cores grow by the accretion of pebble-size objects, as the latter drift 
towards the star due to gas drag. 

Aims. We investigate the dynamics of pebble-size objects in the vicinity of planetary embryos of 1 and 5 Earth masses 
and the resulting accretion rates. 

Methods. We use hydrodynamical simulations, in which the embryo influences the dynamics of the gas and the pebbles 
suffer gas drag according to the local gas density and velocities. 

Results. The pebble dynamics in the vicinity of the planetary embryo is non-trivial, and that it changes significantly 
with the pebble size. Nevertheless, the accretion rate of the embryo that we measure is within an order of magnitude 
of the rate estimated in LJ12 and tends to their value with increasing pebble-size. 

Conclusions. The model by LJ12 has the potential to explain the rapid growth of giant planet cores. The actual accretion 
rates however, depend on the surface density of pebble size objects in the disk, which is unknown to date. 

Key words. Accretion, Planets: formation, hydrodynamics 



1. Introduction 

The formation of the massive cores of giant planets within 
the short timescale allowed by the survival of a proto- 
planetary disk of gas and solids (a few My; Haisch et 
al. 2001) is still an open problem. In the classical view, 
these cores form by collisional coagulation from a disk of 
km-sized planetesimals, through the well-known processes 
of runaway (Greenberg et al. 1978; Wetherill & Stewart 
1989) and oligarchic growth (Ida & Makino 1993; Kokubo 
& Ida 1998). In principle these processes should continue 
until the largest objects achieve an isolation mass, which 
is a substantial fraction of the initial total mass of local 
solids. If the initial disk is sufficiently massive (about 10 
times the so-called Minimal Mass Solar Nebula or MMSN; 
Wcidcnschilling 1977; Hayashi 1981), it is expected that 
cores of <~ 10 Earth masses (Me) form beyond the so- 
called snowline (The orbital radius beyond which tempera- 
ture is cold enough that water condenses into ice; Podolak 
& Zucker 2004), as required in the core-accretion model for 
giant planet formation (Thommes et al. 2003; Goldreich et 
al. 2004; Chambers 2006). 

A-body simulations, though, show that reality is not so 
simple. When the cores achieve a mass of about 1 M® they 
start to scatter the planetesimals away from their neigh- 
borhood, instead of accreting them (Ida & Makino 1993; 
Levison et al. 2010), which slows their accretion rate sig- 
nificantly. It has been proposed that gas drag (Wetherill & 
Stewart 1989) or mutual inelastic collisions (Goldreich et al. 
2004) prevent the dispersal of the planetesimals by damp- 
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ing their orbital eccentricities, but in this case the cores 
open gaps in the planetesimal disk (Levison & Morbidelli 
2007; Levison et al. 2010), like the satellites Pan and Daphis 
open gaps in Saturn's rings. Thus the cores isolate them- 
selves from the disk of solids. This effectively stops their 
growth. It has been argued that planet migration (Alibert 
et al. 2004) or the radial drift of sub-km planetesimals due 
to gas drag (Rafikov 2004) break the isolation of the cores 
from the disk of solids but, again, A-body simulations show 
that the relative drift of planetesimals and cores simply col- 
lects the former in resonances with the latter (Levison et al. 
2010); this prevents the planetesimals from being accreted 
by the cores. In fact, only planetesimals smaller than a few 
tens of meters drift in the disk fast enough to avoid trap- 
ping in any resonance with a growing core (Weidenschilling 
& David 1985). 

In a recent paper, Lambrechts & Johansen (2012), here- 
after LJ12, have proposed a new model of core growth, 
which argues that, if the mass in the disk is predominantly 
carried by pebbles of a few decimeters in size, the largest 
planetesimals accrete pebbles very efficiently, rapidly grow- 
ing to several Earth masses (see also Johansen & Lacerda 
2010; Ormel & Klahr 2010 and Murray-Clay et al. 2011). 
More specifically, this model builds on the recent planetesi- 
mal formation model (Youdin & Goodman 2005; Johansen 
et al. 2006, 2007, 2009) in which large planetesimals (with 
sizes from <~ 100 up to <~ 1,000km) form by the collapse of 
a self-gravitating clump of pebbles, concentrated to high 
densities by disk turbulence and the streaming instability. 
The mechanism by which, once formed, planetesimals can 
keep accreting background pebbles is described hereafter. 
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Pebbles are strongly coupled with the gas; thus they 
encounter the already-formed planetesimals with a velocity 
Av that is equal to the difference between the Keplerian 
velocity and the orbital velocity of the gas (slightly sub- 
Keplerian due to the outward pressure gradient). LJ12 de- 
fine the planetesimal Bondi radius as the distance at which 
the planetesimal exerts a deflection of one radian on a par- 
ticle approaching with a velocity Av: 



Rb = 



GM 



A? 



(1) 



where G is the gravitational constant and M is the planetes- 
imal mass (obviously the deflection is larger if the particle 
passes closer than Rb)- LJ12 showed that all pebbles with a 
stopping time t s smaller than the Bondi time ts — Rb/ Av 
that pass within a distance R = (ts/is) 1 / 2 Rb spiral down 
towards the planetesimal and are accreted by it. Thus, the 
growth rate of the planetesimal is: 



M = irpR 2 Av 



(2) 



where p is the volume density of the pebbles in the disk. 

From (IT]), the Bondi radius grows with the planetesi- 
mal mass. LJ12 also showed that, when the Bondi radius 
exceeds the scale height of the pebble layer, the accretion 
rate becomes 



M = 2i?£Av 



(3) 



where £ is the surface density of the pebbles. Moreover, 
when the Bondi radius exceeds the Hill radius Rh , the ac- 
cretion rate becomes 



M = 2R H Y,v H 



(4) 



where vh is the Hill velocity (i.e. the difference in Keplerian 
velocities between two circular orbits separated by Rh)- 

With these formulae, and assuming that £ stays con- 
stant and is close to the nominal density of solids in the 
MMSN, LJ12 showed that the formation of 10 M E cores is 
possible within 1 My essentially anywhere in the disk (up 
to - 50AU). 

There are two main advantages in the LJ12 model. First, 
it can form 10 Me cores of giant planets within the life- 
time of the disk, a result very difficult to achieve by other 
models. Second, because the accretion rate |2]) is very sensi- 
tive on the planetesimal mass (M oc M 2 ), in practice only 
the largest planetesimals formed in the turbulent model 
can effectively grow in mass by this process: the minimal 
mass for triggering significant Bondi accretion (see eq. ffl is 
about the mass of Ceres in the asteroid belt and about 
the mass of Pluto in the Kuiper belt. Thus this model 
explains the maximal sizes observed in the asteroid and 
Kuiper belt populations. In essence, in this model bodies 
smaller than Ceres (respectively Pluto for the Kuiper belt) 
remained small bodies (the asteroids and KBOs we see to- 
day) , whereas those bigger than this threshold kept accret- 
ing pebbles and became massive objects (embryos) which 
then were removed by migrating away and (possibly) par- 
ticipating to the build-up of the giant planets. Both these 
aspects of the model are very appealing. 

However the study conducted in LJ12, both in the an- 
alytic and in the numerical parts, assumes that the motion 
of the gas is not perturbed by the planetesimal. This as- 
sumption is good for a Ceres-mass planetesimal, accreting 



as in ([2]), but it is far from reality for planetary embryos 
(Earth mass or larger), accreting through their Hill sphere 
as in Q. In fact, these objects modify the gas streamlines 
significantly: a spiral density wave is formed in the disk and 
the gas near the orbit of the planet has horseshoe motion. 
An over-density of gas is also formed inside the planet's 
Hill sphere. It is not clear a priori what are the effects of 
these structures on the pebble accretion rate. This is pre- 
cisely what we investigate in this paper with more realistic 
hydro-dynamical simulations. 

In section 2, we describe our methods: the simulation 
tool that we have developed and the parameters that we 
adopt. In section 3 we present our results, for two embryo 
masses and 4 pebble sizes. Our goal is three-fold: (i) de- 
scribe and understand the dynamics of the pebbles for the 
different mass and size cases; (ii) evaluate the accretion rate 
by the embryo and compare it with the LJ12 estimate and 
(iii) evaluate the "filtering factor" , that is the fraction of the 
pebbles that do not drift by the orbit of the planet because 
they are accreted by the embryo. This factor is important. 
If it is large, of a sequence of embryos radially distributed 
in the disk, only the outermost one(s) can accrete; instead, 
if it is small then the full system of embryos can grow, in 
an oligarchic fashion. Our conclusions and discussion of a 
coherent scenario of giant planet formation conclude the 
paper in section 4. 

2. Methods 

Our simulation software is based on the hydro-dynamical 
code FARGO, developed to study planet-disk inter- 
actio ns in Masset (2000) and p ublicly available at 
http://fargo.in2p3.fr/spip.php7auteurl. In that code, how- 
ever, the N-body interaction among the bodies in the sys- 
tem is studied with a Runge-Kutta integrator of 5th order. 
Because the time-step of the integration is fixed, the Rungc- 
Kutta algorithm is not adequate to treat accurately the 
close encounters between objects (e.g. the encounters be- 
tween the pebbles and the planetary embryo). Also, there 
is no prescription in the original code to detect collisions. 

Thus, we replaced the Runge-Kutta subroutine in 
FARGO with a code known as Symba. The latter is a 
variable-timestep symplectic code developed in Duncan et 
al. (1998) to simulate quasi-Keplerian N-body systems with 
mutual close encounters. Symba also identifies collisions 
and merges the bodies that collide. A technical difficulty 
in interfacing FARGO with Symba was that the former is 
written in C-language and the latter in Fortran. This re- 
quired extensive modifications of several subroutines of the 
FARGO package. 

In our simulations, the embryo is a massive object: it 
influences the evolution of the gas but, for simplicity, we 
cancel the influence of the gas disk on the embryo, so that 
the latter remains on a fixed, non-migrating orbit. This 
approximation is reasonable, as long as the migration of the 
embryo is slow compared to the radial drift of the pebbles 
due to gas drag. 

In the FARGO code, the gas is discretized over a po- 
lar two-dimensional grid. However, in our modified code, 
the pebbles can evolve in a three-dimensional space. The 
gas drag on the pebbles is then computed as follows. We 
consider a spherical coordinate system r, 9, <f>, where r, 6 
are the radial and azimuthal coordinates on the disk mid- 
plane. At each timestep we identify the disk grid cell where 
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the spherical projection of the pebble falls on the midplanc. 
This depends on the r, coordinates of the pebble. We con- 
sider the values of the gas surface density E and radial and 
tangential velocities v r ,vg that the gas has in that cell as 
well as in its 8 neighboring cells. We interpolate the set of 
9 values for each field with a polynomial function of typej 
(a + bx + cx 2 ) + (a' + b'x + c'x 2 )y + (a" + b"x + c"x 2 )y^ 
and we use it to evaluate E, v r and vg of the gas at the r, (t 
location of the pebble. Moreover, we assume that the verti- 
cal velocity of the gas v z is zero. This sets the local relative 
velocity 5v of the gas and the pebble. The local volume 
density of the gas at the pebble location is computed as 
p = E/(v / 7ri7) cxp(— z 2 /H 2 ) where z is the vertical coordi- 
nate of the pebble and H/r is the assumed aspect ratio of 
the disk, given as an input of the simulation. Finally, the 
drag suffered by the pebble is computed as v = —lSv/t s , 
where t s is the stopping time. The stopping time is also 
evaluated locally following Adachi et al. (1976), given the 
size and bulk density of the pebble (given as an input of 
the simulation) and the local volume density of the gas and 
Mach & Knudsen numbers. 

2.1. Simulation parameters 

Our simulations follow the dynamics of one embryo and a 
set of pebbles. We consider planetary embryos of 1 or 5 
Me- They are placed on a circular non-migrating orbit at 
0.8 AU. The disk's surface density is equal to 1800g/cm 2 
at 1 AU, corresponding to the MMSN, with a radial profile 
proportional to 1/r. The aspect ratio is 5.6%. We follow an 
a-prescription for the disk's viscosity (Shakura & Sunyaev 
1973), with a — 6 x 10 -3 . The parameters above consti- 
tute just one set of (typical) disk parameters. Different disk 
parameters would affect mainly the pebble stopping time. 
Given that we will study the dynamics of pebbles over a 
wide range of sizes, our exploration of the stopping-time 
parameter space will be exhaustive even if working with 
just one disk model. 

The disk is modeled over a grid extended from 0.75 to 
1.12 AU, with a resolution of 160 concentric rings and 720 
sectors. As the disk radial extension is narrow, we use non- 
reflecting boundary conditions, so that the spiral density 
wave launched by the embryo does not make the surface 
density of the disk rough with multiple reflections. 

An image of the surface density of the disk is shown in 
Fig[l] with a color scale. The black curves show stream- 
lines of the gas in the reference frame corotating with the 
embryo. The position of the embryo is highlighted by a lo- 
cal maximum of the disk's surface density. Notice that the 
spiral density wave (the slanted over-density structure de- 
parting from the planet) does not have any reflection at 
the border of the frame. Therefore, it is not necessary to 
simulate a wider radial portion of the disk. 

The gravitational interaction between the embryo and 
the gas is smoothed as 1/(A + e) 2 where A is the distance 
between the embryo and a gas fluid element and e is a 
smoothing parameter. We adopt the quite standard choice 
e = 0.6R H . 

We simulated pebbles of 5cm, 20cm, lm and 10m in ra- 
dius, with a bulk density of lg/cm 3 . With our disk model 
this corresponds to the following stopping times at 1 AU: 
0.017, 0.22, 2.78, 74.07 cu respectively, where cu is the time 
code unit, which is 1/27T of an year. Even if some of these 
particles are rather boulder-sized, we still call them peb- 
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Fig. 1. Color map of the surface density of the gas in a disk in 
the presence of a 5Me embryo at r — 0.8, 9 = 3.22. The black 
curves show some gas streamlines and the blue, red, yellow and 
white curves depict the trajectories of 5cm, 20cm, lm and 10m 
pebbles, respectively, in the frame corotating with the embryo. 

bles, for simplicity. The pebbles are assumed to orbit on 
the same plane as the disk and the embryo so that our sim- 
ulations become effectively two-dimensional. A trajectory 
of one pebble of each size is shown in Fig. [I] and described 
in the next section. The timestep used in the integration is 
the minimum between 1 /4 of the particle stopping time and 
0.1 cu. The latter is the timestep dt imposed by the CFL 
condition (Courant et al. 1928) for the numerical solution 
of the hydrodynamical equations, given the resolution and 
the extension of our grid. Remember that the Symba algo- 
rithm effectively reduces the timestep for the integration of 
the particles approaching the embryo, through a clever sub- 
division of the embryo's gravitational potential (Duncan et 
al. 1998). 

To measure the accretion rate and the filtering factor 
(see last paragraph of Sect. 1) we proceed in two steps. We 
first integrate one pebble starting from the outer boundary 
of the disk at r = 1.12 AU at opposition with respect to 
the embryo. As the pebble migrates inward, we record the 
sequence of heliocentric distance values that the pebble 
has at subsequent oppositions. We denote by the small- 
est of the r-fc values with > 0.8 AU, i.e. the heliocentric 
distance at the last opposition before crossing the orbit of 
the embryo^ Then, we simulate a large number of pebbles 
(usually 100, but this number can change from case to case 
to achieve an appropriate resolution), initially at opposi- 
tion and with heliocentric distances uniformly distributed 
in the range [rjy, rjy-i]- Of this beam of pebbles, we mea- 
sure the fraction F that are accreted by the planet. F is the 
"filtering factor" defined in the introduction, while 1 — F is 
the fraction of the beam that drift across the orbit of the 
embryo. The accretion rate is then: 

M = E(rjv-i - r N ) x FAv (5) 

where Aw is the difference in orbital velocity between the 
embryo and a pebble placed in the middle of the [rjv, FjV— l] 
interval and E is the surface density of the disk of pebbles. 

1 For the 10m "pebble", we chose instead N = 2, so that it 
is still far enough from the embryo to preserve a quasi-circular 
orbit. 
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Thus, our result is parametric in E and can be directly 
compared with the expression Q from LJ12. 

3. Results 

Fig. [T] shows the trajectory of one pebble for each of the 
4 sizes that we consider in this study. The evolution is 
shown in a rotating reference frame, with the embryo fixed 
at r = 0.8 AU and 9 = 3.22. Only the last synodic revolu- 
tion is shown, i.e. that leading to the crossing of the orbit 
of the embryo. The trajectories of the 5cm, 20cm. lm and 
10m pebbles are depicted in blue, red, yellow and white, 
respectively. The mass of the embryo is 5 Me. 

Several considerations should be made, by looking at 
these trajectories. First, notice that the stopping time 
(which increases with particle size) is not simply related 
to the migration speed. The particles with the fastest ra- 
dial migrations are those with radii of 20cm and lm (with 
migration rates roughly equal to each other). This is well 
known (Weidenschilling 1977b). The 5cm particle migrates 
more slowly because it is more coupled with the gas (which 
does not have any net radial motion); the 10m particle mi- 
grates more slowly than the lm particles because it is less 
sensitive to gas drag. 

The U-turns drawn by the trajectories of the 20cm and 
lm pebbles are not related to horseshoe dynamics. They 
are simply due to the fact that the angular velocity of a 
Keplerian orbit decreases as 1/r 3 / 2 so that, in a reference 
system corotating with an object at a fixed distance r (like 
Fig.[l]corotating with the embryo at ro = 0.8AU) a particle 
moves from the right to the left if it has r > ro and from 
the left to the right if it has r < r . Thus, as the particle 
drifts from r > ro to r < ro it has to make a U-turn. In 
other terms, these bended trajectories are not due to the 
gravitational effect of the embryo but just to the rigid rota- 
tion of the reference frame. Only particles passing close to 
the embryo (shown below) are affected by the gravitational 
attraction of the latter. 

Instead, the U-turn of the 5cm pebble is due to horse- 
shoe dynamics. This is apparent from the fact that the ra- 
dial motion is fast only in the vicinity of the embryo. The 
particle shown in Fig. [I] traces approximately the maximal 
width of a particle horseshoe trajectory. Notice that the 
width of the horseshoe region for the particles is much nar- 
rower than that for the gas, highlighted by the U-turning 
streamlines. 

Finally, the trajectory of the 10m object shows radial os- 
cillations as it approaches towards the embryo. This is due 
to its orbital eccentricity, acquired during previous passages 
at conjunction with the embryo, and not fully damped by 
the gas drag in half a synodic period (i.e. from conjunction 
to opposition) . This is not the case from the other particles 
for two reasons: gas drag is stronger and erases the eccen- 
tricity faster (e.g. the 5cm pebbles) and/or given their fast 
radial drift the pebbles passed too far from the embryo at 
the previous conjunction to acquire a large eccentricity (e.g. 
the 20cm and lm particles). 

We now focus on the dynamics in the very vicinity of 
the embryo, leading to pebble accretion and the embryo's 
growth. 

Fig. [2] shows the trajectories of several 5cm pebbles as 
they pass close to the IMe embryo. Each line is a differ- 
ent particle. The x-axis reports the difference 59 between 
the 9 values of the pebble and the embryo. Notice on the 
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Fig. 2. the trajectories of 5cm pebbles approaching a 1 Me 
embryo in relative r, 6 coordinates. 
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Fig. 3. The same as Fig. |2| but for 20cm pebbles. 



right hand side of the figure the particles doing the horse- 
shoe U-turn (those which always have 59 > 0). Instead, 
the particles that are initially in circulation, i.e. that have 
r > 0.8 AU and 59 passing from positive to negative, even- 
tually cross the embryo's orbit far behind the embryo it- 
self, when 59 < —0.6. Notice that the U-turn of the pebbles 
occurs at r ~ 0.798, i..e slightly inside the embryo's or- 
bit (r = 0.8). This is because these pebbles are strongly 
coupled with the gas, which has a sub-Keplerian orbital ve- 
locity; thus the exact corotation radius is shifted towards 
the Sun. In total, 7 pebbles collide with the embryo, as- 
sumed to have the Earth's radius, as one can see from 
the trajectories diving towards the point with coordinate 
(0, 0.8). Six of these pebbles collide with the embryo in the 
approach phase, when 59 is evolving from positive to nega- 
tive. However, one pebble hits the embryo "from the back" , 
i.e. after having crossed the embryo's orbit and with 59 
evolving from negative to positive. 

Fig. [3] is the same but for 20cm pebbles. The figure is 
qualitatively similar to the previous one, but the trajec- 
tories are much more inclined due to the fact that these 
pebbles drift much faster towards the Sun (beware that 
the scales in Fig. [2] and [3] are different) and, consequently, 
the particles that cross the orbit of the embryo at negative 
59 do so closer to the embryo (with 59 > —0.2 instead of 
< —0.6). The thick ellipse drawn around (0,0.8) marks a 
circle centered on the embryo with a radius equal to 1/3 
of its Hill radius. The resolution of the output of the sim- 
ulation is too coarse to resolve the trajectories inside this 
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Fig. 4. The same as Fig. [2] but for lm particles. 




Fig. 5. The same as Fig. [2] but for 10m particles. 



radius. However, the internal timestep of the Symba algo- 
rithm is good enough to resolve the physical collision of the 
pebbles with the embryo. There are 18 pebbles entering in- 
side the circle, 4 of which enter "from the back" , i.e. after 
having crossed the embryo's orbit. Of these 18 pebbles, 13 
physically collide with the embryo. The remaining ones are 
slingshot away by the embryo during the close encounter 
and exit the domain covered by our grid. 

Fig. [4] is for one-meter particles. The major difference 
with respect to the previous figure is the strong kink that 
particles receive when passing from positive to negative 59. 
This happens because these particles are less coupled to 
the gas and therefore they can be scattered to large ec- 
centricity when they pass in conjunction with the embryo. 
Consequently, the 1-m particles pass much further away 
from the embryo than the 20cm-particles when they cross 
the embryo's orbit, and therefore they cannot collide with 
the embryo "from the back" . 

Fig.[5]is for lOm-particles. For these particles we change 
the scale of the plot, showing the whole disk up to r = 
0.86 AU. This is because gas drag is weak, and the particles 
conserve part of the eccentricity that they acquire during 
the previous conjunction with the embryo for a full syn- 
odic period. Consequently, their dynamical behavior when 
crossing the embryo's orbit inherits the dynamical history 
recorded over several synodic revolutions. Because the ra- 
dial drift during a synodic revolution is small, all particles 
have to pass eventually very close to the embryo, so that 
many of them (40%) are accreted. The particles that safely 
cross the embryo's orbit do so with a horseshoe U-turn at 
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Fig. 7. The accretion rate as a function of pebble size for em- 
bryos of IMe (filled dots) and 5Me (open dots) respectively. 
The dotted and dashed lines indicate the accretion rates esti- 
mated in LJ12 for these two embryo masses. 



59 > 0. If the embryo is more massive the fraction of ac- 
creted particles increases and, for a 5 Me embryo, the ac- 
cretion efficiency is 100% (see Fig. |6j. 

We now come to a synthesis of the results for what con- 
cerns the embryo's mass accretion rate and filtering factor. 

Fig. [7] shows the mass accretion rate for the embryos of 
masses IMe (filled dots) and 5Me (open dots) respectively, 
as a function of the pebbles size. The mass accretion rate is 
normalized by £ (the surface density of the pebbles). The 
dotted and dashed horizontal lines show the accretion rate 
estimated by LJ12, for embryos of masses IMe and 5 Me 
respectively. As one can see, our numerically measured ac- 
cretion rates are within an order of magnitude of the es- 
timated rates, and tend asymptotically to the estimated 
values for growing pebble sizes. The drop of the accretion 
rate with particle size is due mainly because a larger frac- 
tion of the particles that enter the Hill sphere of the embryo 
are dragged away by the flow of the gas, instead of falling 
onto the central object. This phenomenon was expected in 
LJ12, who predicted that their accretion rate is valid for 
pebbles with stopping time of the order 0.1-1. Remember 
that 20cm pebbles in our simulation have t s = 0.22. Thus, 
our results suggest that the nominal accretion rate of LJ12 
applies instead from t s ~ 1; however, it remains valid up 
to stopping times larger than predicted, i.e. up to at least 
~ 100 (the value of t s for the 10m particles is 75). 
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Fig. 8. The filtering factor as a function of pebble size for em- 
bryos of 1Mb (filled dots) and 5Me (open dots) respectively. 



The scaling relative to the embryo's mass predicted by 
LJ12 (M oc RhVh, i-e. M oc M 2 / 3 ) is confirmed by our 
numerical simulations for all tested pebble sizes. 

Figure [8] shows the filtering factor as a function of pebble 
size, again for 5Me (open dots) and lMg (filled dots) em- 
bryos. The filling factor scales again as Af 2 / 3 . The filtering 
factor decreases with increasing drift speed and therefore it 
is the smallest for the m-sized particles. Except for parti- 
cles of multiple meters, the filtering factor is small (a few 
percent to 10%), which implies that a system of numerous 
embryos (with a number of objects roughly proportional to 
1/F) can grow simultaneously in the disk. Notice that the 
embryos should grow in an oligarchic fashion, due to the 
M 2 / 3 dependence of the accretion rate. Instead, if the disk 
is dominated by multi-meter particles, only a few embryos 
can grow, with eventually the outermost one capturing all 
the material. 

In the LJ12 model, there is no upper limit in mass for 
an accreting embryo, as long as pebbles are available: an 
embryo keeps accreting at a rate proportional to Af 2 / 3 £. 
In reality, however, if the total mass of the planet (the mass 
of the solid embryo plus the mass of the gas eventually ac- 
creted in its atmosphere) becomes large enough, the surface 
density of the gas distribution is modified. The neighbor- 
hood of the planet orbit becomes partially depleted of gas 
(a shallow gap is opened) and this changes the pressure 
gradient in the disk and the angular velocity of the gas. 
Eventually the gas becomes super-Keplerian at some loca- 
tion outside of the planet's orbit and this stops the inward 
radial drift of pebbles and small planetesimals. The accre- 
tion of the planet by the LJ12 mechanism suddenly stops. 
Fig [9] shows the ratio between the azimuthal velocity of the 
gas and that of a Keplerian circular orbit as a function 
of radius, for planet masses of 30, 50 and 100 Me placed 
at r — 0.8. For these tests, we have increased the disk's 
radial extension to the interval 0.5-1.3 AU. For the disk 
parameters that we chose (see sect. 2.1), we find that the 
disk becomes locally super-Keplerian is the planet's mass 
is ~ 50Me- This mass would be reduced if the disk had 
a lower viscosity and/or scale height. Whatever the exact 
value of the limit mass of the planet, our result shows that 
the LJ12 mechanism is capable of feeding a solids to a giant 
planet until its mass is several tens of Earth masses. 

This result has two main implications. First, the run- 
away accretion of the atmosphere onto the embryo might 
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Fig. 9. The gas azimuthal velocity relative to the Keplerian 
value, as a function of orbital radius, for disks with embed- 
ded planets of 30, 50 and 100 Earth masses. The horizontal 
line marks the boundary between sub-Keplerian and super- 
Keplerian rotation. 



start earlier than in the classic Pollack et al. (1996) model. 
In fact, in Pollack et al. the embryo's accretion stalls when 
the mass of the embryo is about 10 Me, because the feed- 
ing zone for solids is strongly depleted. Then, the embryo's 
accretion of both gas and solids continues very slowly, with 
the consequence that the critical mass for triggering run- 
away capture of gas (~ 30Mb) is reached only after sev- 
eral My. In the LJ12 scenario, the mass accretion of the 
embryo would not slow down, so that the runaway accre- 
tion of gas, in principle, might start at an earlier time. 
We notice however that a high accretion rate of solids has 
also a negative effect on gas accretion because it releases 
a lot of energy on the embryo which needs to be evacu- 
ated (Dodson-Robinson & Bodenheimer 2010); thus, the 
net effect of pebble-accretion on gas-accretion will need to 
be investigated in hydrodynamical simulations accounting 
for heat transfer. Second, the accretion of pebbles until the 
planet mass achieves ~ 50 Me can help explaining the large 
ratio between heavy elements and hydrogen and helium in 
the atmosphere of Jupiter, which is 3-4 times solar (Wong 
et al. 2004). Enriching the giant planets in heavy elements 
is not straightforward (Guillot & Gladman 2000) by other 
mechanisms; thus this is another appealing aspect of the 
LJ12 model. 



4. Conclusions 

In this paper we have tested the scenario of giant planet 
core formation proposed by LJ12 with hydrodynamical sim- 
ulations that fully account for (i) the interaction between 
the growing core and the gas of the disk and (ii) the local 
drag exerted by the gas on pebbles and boulders. 

We have found that the pebble dynamics in the vicinity 
of the planetary embryo is non-trivial, and that it changes 
significantly with the pebble size. Nevertheless, the accre- 
tion rate of the embryo that we measure is within an order 
of magnitude of the rate estimated in LJ12 and tends to 
their value with increasing pebble-size. 

The accretion of pebbles can continue until the embryo's 
mass is of the order of 50 Earth masses (in solids and gas) . 
This may have important implications on the onset of run- 
away accretion of gas by the growing giant planets and 
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can help explaining the enrichment in solids observed in 
Jupiter's atmosphere. 

The actual accretion rate of an embryo depends on the 
amount of mass E available in pebbles, which is not known 
a priori, given that pebbles are consumed by the forma- 
tion of planetesimals and by the accretion of the embryos 
themselves. Nevertheless, the accretion rates that we find 
in Fig. [7] are potentially large. For instance, if a MMSN of 
solids (20g/cm 2 at 1 AU) were available in 20cm pebbles, 
the mass doubling time for a 1 Me embryo would be only 
5,500 years! This illustrates the importance of the LJ12 
model for the growth of giant planets cores. 

Given that the LJ12 model is built in the same frame- 
work as the model that explains the rapid formation of 
planetesimals (Johansen et al. 2007), we believe that the 
community has now, for the first time, a coherent scenario 
to explain the the early phases of planet growth. However, 
we do not see any evident reason for which only a small 
number (4-6) of giant planet cores should grow in the disk, 
as suggested by the number of giant planets of our so- 
lar system, including rogue ice-giants potentially lost dur- 
ing the dynamical evolution that followed planet formation 
(Nesvorny 2011; Nesvorny & Morbidelli 2012). Thus, we 
think that, most likely, the LJ12 model explains the forma- 
tion of massive planetary embryos of a few Earth masses, 
but an additional stage is needed to form the giant planet 
cores (10 Earth masses or more). 

This is the scenario that we envision. The embryos 
formed by the LJ12 mechanism, once massive enough, start 
to migrate in the disk due to planet-disk interactions. 
Recent results on migration in radiatively cooling disks 
(Paardekooper & Mallema 2006; Baruteau & Masset 2008; 
Kley & Crida 2008; Paardekooper et al. 2010; Masset & 
Casoli 2010; Bitsch & Kley 2011) show that the embryos 
migrate from all directions toward an orbital radius where 
migration is canceled out by the compensation of compet- 
ing torques. This convergent migration towards the same 
region can promote the mutual accretion of embryos, even- 
tually reducing a system of a large number of embryos into 
a system of a smaller number of larger objects, i.e. a hand- 
ful of giant planet cores (Horn et al. 2012). Admittedly, 
this scenario is still speculative and more work is needed 
to prove its validity. We stress, however, that a final phase 
of core formation characterized by mutual collisions of em- 
bryos would explain, in a natural way, the massive impacts 
that are needed to explain the current obliquities of Uranus 
and Neptune (Morbidelli et al. 2012). 
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